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Quantum phase estimation algorithm has been successfully adapted as a sub frame of many other 
algorithms applied to a wide variety of applications in different fields. However, the requirement of 
a good approximate eigenvector given as an input to the algorithm hinders the application of the 
algorithm to the problems where we do not have any prior knowledge about the eigenvector. 

In this paper, we show that the principal eigenvalue of an irreducible non-negative operator can 
be determined by using an equal superposition initial state in the phase estimation algorithm. This 
removes the necessity of the existence of an initial good approximate eigenvector. Moreover, we 
show that the success probability of the algorithm is related to the closeness of the operator to a 
stochastic matrix. Therefore, we draw an estimate for the success probability by using the variance 
of the column sums of the operator. This provides a priori information which can be used to know 
the success probability of the algorithm beforehand for the non-negative matrices and apply the 
algorithm only in cases when the estimated probability reasonably high. Finally, we discuss the 
possible applications and show the results for random symmetric matrices and 3-local Hamiltonians 
with non-negative off-diagonal elements. 


I. INTRODUCTION 

In the recent decades, there have been many quantum algorithms proposed for the problems inefficient to solve on 
classical computers. The quantum phase estimation algorithm [1, 2] is a leading illustration of such algorithms. It is 
used for quantum simulations to estimate the eigenenergy corresponding to a given approximate eigenvector of the 
unitary evolution operator of a quantum Hamiltonian. Furthermore, with different settings, it has been adapted as 
a sub frame of many quantum algorithms applied to wide variety of applications in different fields (see the review 
article ref.[3] and the references therein). However, the success of the algorithm partly depends on the pre-existing 
approximate eigenvector. Thus, it may fail to output the right eigenvalue in the cases where the approximation to 
the corresponding eigenvector is not good enough. This hinders the uses of the algorithm when there is not enough 
prior knowledge to procure a good approximation to the eigenvector. 

A symmetric matrix is called non-negative if all of its elements are greater than or equal to zero. It is irreducible 
if the corresponding graph is strongly connected or if it is not reducible to a block-diagonal form by any column-row 
permutation. Due to Perron-Frobenius theorem [4], an irreducible non-negative matrix have a positive eigenvector 
(all elements are positive) with an associated positive principal eigenvalue whose magnitude is greater than the rest 
of the eigenvalues. These matrices have been studied extensively [5, 6], the distribution of the coefficients of their 
principal eigenvectors have been related to the matrix elements [7], and the sum of the coefficients have been shown 
to be related to the number of walks in molecular graphs [8]. 

The phase estimation algorithm (PEA) mainly uses two quantum registers:jregl) initially set to zero state and \reg 2 ) 
holding an initial approximate eigenvector. In this article, we show that for irreducible non-negative matrices, one can 
obtain the principal eigenvalue by using an equal superposition input state in PEA: i.e., instead of an approximate 
eigenvector, initial value of \reg 2 ) is set to | 0 ) and then put into equal superposition state by applying Hadamard 
gates. In the output of the algorithm, this generates each one of the eigenvalues with the probability determined 
by the normalized sum of the coefficients of the associated eigenvector. In addition, because all eigenvectors but the 
principal one include positive-negative elements, we show that in most of the random cases the probability to see 
principal eigenvalue in the output is much larger than the others. In also some cases, we show that first applying 
Hadamard gates to the second register in the output (or measuring second register in the Hadamard basis) and then 
measuring the first register when the measurement outcome of the second register is in |0.. .0) state increase the 
success probability further in the output. 

Eigenvector-coefficients of stochastic matrices sum to zero for all but the principal eigenvector. Therefore, for these 
matrices one can generate the eigenvalue in the phase estimation algorithm with probability one [9]. Since any given 
symmetric irreducible non-negative matrix can be converted into a stochastic matrix by a diagonal scaling matrix; 
using the closeness of the matrix to a stochastic matrix, we also show that the success probability of the algorithm can 
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be predicted beforehand. We finally compare the estimated success probabilities with the computed ones for random 
symmetric matrices and 3 -local Hamiltonians of different dimensions with non-negative off-diagonal elements. 

In the following sections, we shall first discuss PEA in detailed steps, then draw the estimates for the success 
probability and finally show the possible applications of the algorithm. 


II. QUANTUM PHASE ESTIMATION ALGORITHM 

The phase estimation algorithm (PEA) in general sense finds the value of (pj for a given approximate eigenvector 
l/ij) in the eigenvalue equation C/|/ij) = The algorithm mainly uses two quantum registers: Viz. |regl) 

initially set to zero state and \reg2) holding an approximate eigenstate of a unitary matrix U. The first operation in 
the algorithm puts \regl) into the equal superposition state. In this setting, a sequence of operators, U'^\ controlled 
by the jth qubit of \regl) are applied to \reg2)- Here, j = 1 .. .m and m is the number of qubits in \regl) and also 
determines the precision of the output. These sequential operations generate the quantum Eourier transform (QFT) 
of the phase on \regl). Therefore, the application of the inverse quantum Fourier transform {QFT'^) turns the value 
of \regl) into the binary value of the phase. Consequently, one measures \regl) to obtain the phase. Here, if the 
unitary operator U is the time evolution operator of a quantum Hamiltonian TL, i.e. U = then one also obtains 
the eigenenergy of that Hamiltonian. 

For a symmetric irreducible non-negative matrix "H of order 2 ” with ordered eigenvalues as > • • • > 02" and 
associated eigenvectors \gi)... \g,2"), it is known that |/ri) is the only eigenvector with all positive coefficients. Assume 
unitary operator U = and its powers with j = 0 .. .m are readily available for PEA. To estimate the value of 
01 on the first register, in our setting, we simply modify the conventional phase estimation algorithm in the following 
way: 

• Instead of an approximate eigenvector, initial value of \reg2) is set to | 0 ) and then put into equal superposition 
state by applying Hadamard gates. 

• (optional) In the output, we change the basis of the second register by applying the Hadamard gates. Then, 
if the measurement of the second register is equal to |0) state, then the principal eigenvalue is estimated on the 
first register. This modification increases the success probability of the algorithm further when the probability 
of measuring the principal eigenvalue is higher then the other eigenvalues. 

In the following subsection, we shall describe the phase estimation algorithm in steps, also shown in Fig.l, to show 
how the above modifications have an impact on the algorithm. 


A. Steps of the Algorithm 

Here, the algorithm is described in details by showing how the quantum state changes in each step: 

• The system is prepared in a way so that it includes two register: viz., \regi) and \reg2) with m and n number 
of qubits, respectively. 

• Both registers are initialized into zero state: \^o)=\regi)\reg 2 )=\ 0 )\^)■ 

• Then, the Hadamard operators are applied to both registers: 

2™+"-l 

|0i) = (i7®™®i7®")|O)|O) = -^ |x) (1) 

• As in the customary phase estimation algorithm; the unitary evolution operators C/^'* controlled by the jth qubit 
of \regi) are applied to \reg2), and finally the inverse of the quantum Fourier transform (QATl) is applied to 
\regl). 

At this point of the algorithm, we have a quantum state in which \regi) and \reg2) hold the superposition of 
the eigenvalues and the associated eigenvectors, respectively: 


N-l 

\ip2) = ^ aj |0j) \gj). 
j=o 


( 2 ) 
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Here, the amplitude aj is related to the angle between the equal superposition state and the eigenvector |/r_,) and 
can be described as the normalized sum of the coefficients of the eigenvectors since \reg 2 ) was 1^) 

the beginning: 


2 ” 

aj = (3) 

i=l 

If \reg 2 ) is written in the Hadamard basis, |H-+ +) + |H-H—)H-h |--), where |±) = l/-\/(2)(|0) ± 

|1)); then the following quantum state is obtained: 


IV' 2 ) — ^ /3ii \4>j) |h-1") 


i=i 


H-h ^ /32"i I-) 

i=i 


(4) 


where j3ijS are new coefficients. Now, since there is only one eigenvector, viz. |/ii), with all positive real elements, 
in the Hadamard basis, we can expect this positive eigenvector to be the closest state to | + + • • • + +). 

• To revert the above state back to the standard basis, the Hadamard operator iJ®” is again applied to \reg 2 )'- 

2” 2" /g\ 

= ^ Pij \4>j) |0 ... 0) H-h ^ (32^j \(t)j) |1 ■ • • 1) 

i=i i=i 


• \reg 2 ) is measured in the standard basis. As a result, for \reg 2 )= |0.. .00), the system collapses to the state 
where the phase associated with the positive eigenvector is expected to be highly dominant in the first register. 

• In the final step, the first register is measured to obtain the phase (f)i and get the eigenvalue of %. These steps 
also drawn in Fig.l. 



FIG. 1. Circuit design for the modified phase estimation algorithm. In the end, note that the application of the Hadamard 
gates to the second register and the measurement on this register is optional. 


B. The Success Probability 

It is easy to see that the success probability of the algorithm without the optional Hadamard gates on \reg 2 ) in 
the output is determined from the normalized sum of the coefficients of the eigenvectors given in Eq.(3). Therefore, 
the success probability for the dominant eigenvalue is af (note that the coefficients of |^i) are all positive, hence, 
a? = |aip.). 

With the Hadamard gates in the output, we first find the probability to measure \reg 2 ) in | 0 ) state and then find 
the probability to see the dominant eigenvalue in \regi) knowing \reg 2 ) = |0): As obvious in Eq.(5), the success 

probability of getting | 0 ) is addition, after measuring | 0 ) on \reg 2 ), the probability to measure 
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\(j)i) on \regi) is |/3iiP- Since we apply an equal superposition as an input, the component of an eigenvector in this 
direction determines the probability to get the corresponding eigenvalue on \regi). More formally, after the application 
of QFT', we get IV' 2 ) = 'l2j=o ^3 I'^i) l^i)’ where aj = l/v2” X]i=i ^ 652 ) of the state in Eq.(4) is measured in 

the Hadamard basis, the probability to see |H-h) is the sum of the amount of the components of the eigenvectors 

in the direction of the initial vector: 


Preg2 — ^ 
i=i 


( 6 ) 


This is also equal to i-O- the probability of measuring \reg 2 ) in |0...0) state at the end. Preg 2 takes the 

smallest value only when all |q;j|s are equal to l/-\/2”. Moreover, when \reg 2 ) = \^ ■ ■ -0), the probability to see on 
\regi) is: 

Preg,=\Pll\'^ = ^- (7) 

^reg2 

Although we have drawn the equalities for probabilities, without knowing the eigenvectors, it is not possible to 
compute exact la^js and so Pregi and Preg 2 - Therefore, we shall try to estimate them: Since ai is the normalized sum 
of the vector elements, it is easy to see that ai > where the equality occurs only when the magnitude of the 

eigenvectors are the columns of the identity matrix. A similar observation is also made in ref. [9] where the principal 
eigenvector is found for a given eigenvalue equal to 1. 

Here, we will relate the matrix "H to a stochastic matrix in order to develop some intuition for the estimation of the 
probabilities. It is known that stochastic matrices have only one eigenvector with coefficients summing to a nonzero 
value. A matrix can be made stochastic by a column wise scaling T-UD, where 2? is a diagonal matrix whose diagonal 
elements are the inverses of the sums of the columns of T-L. The closeness of the matrix Ti, to T-LV may provide a 
prediction to estimate the success behavior of the phase estimation algorithm with an initial superposition state. The 
closeness of T-L to TiV can be defined by \\Ti ~ TiT>\\ where || • || is any matrix norm. This also defines a perturbation 
error. If we normalize this error term by 112211, we get the relative error: 


ei 


1122-222211 

II72II 


< \\p-n, 


( 8 ) 


where X is an identity matrix. We can also look at the relative error in the inverses: 


ei 


1122-1 - 25 - 122-1 

l|72-i|| 


< 11^-77-^11. 


(9) 


Since the matrix 22 also changes the eigenvector elements; when the variance of the diagonal elements of 22 is small, 
we can expect the eigenvectors of 22 to show the similar behaviors as of 2222 and have one eigenvector whose elements 
sum to a much greater value than the others. When 22 = 2222 or 22 = I, left and right principal eigenvectors of 
symmetric 22 are the same and have the elements equal to in the normalized space. Therefore, the variances ui 
and (72 of the diagonal elements of 22“! and 22, the column sums of 22 and their inverses, can somehow describe how 
much the elements of the principal eigenvector deviate from which is also the expectation value for randomly 
generated normalized uniform N numbers. Using these intuitions, we define the estimate Preg 2 Eis: 


p — 

J rPHn - 


Ai + A 2 


( 10 ) 


with 


Ai 


jy_^ 


A 2 


JV_^ 

^ +CT 2 ’ 


( 11 ) 


In addition, af = 1 for stochastic matrices and so for 22 when 22 = 2222. Since the any deviation of the elements of 
the principal eigenvector from would change the equality af = 1, we multiply variance by JV to get total deviation 

from 1. Therefore, we define the estimate af as: 


,2 _ (I — JVai) + (I — 2VCT2) 


( 12 ) 
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Note that af and Preg 2 are only predictions and may fail to provide good estimates for af and Preg 2 j respectively, 
for particularly some structured matrices and matrices whose eigenvectors are close to the columns of an identity 
matrix. For instance, consider the following matrix: 


/21.8214 

0 

0.6118 

0.4983' 

0 

14.2944 

0.4983 

0.6118 

0.6118 

0.4983 

12.1626 

0 

\ 0.4983 

0.6118 

0 

5.4111, 


(13) 


The eigenvalues of this matrix are 5.3537,12.0193,14.4411, and 21.8753, respectively, associated with the following 
eigenvectors: 


/ 0.0304613 
0.0686662 
-0.0077623 
V-0.9971443 


0.0597207 

0.2074209 

-0.9761393 

0.0237068 


0.0215934 0.9975166\ 
-0.9758165 0.0066086 I 
-0.2076080 0.0631720 
-0.0649217 0.0304360/ 


(14) 


As computed from the above, the magnitudes of the sums of the components of the eigenvectors are very close to 
each other: | — 0.90578|, | — 0.68529|, | — 1.22675|, and |1.09773|, respectively. In addition, the second eigenvector, not 
associated to the principal eigenvalue, has the largest sum. 

Also note that one can also predict af from the distribution of the matrix elements. For instance [7, 10], 


maxij 


Mh ^ maxijhjj 


tniTti jh'ij 


(15) 


where hijS are the nonzero matrix elements of %. However, these bounds are generally not tight enough to give good 
estimates in many cases. 

In Sec.IV, we shall show the comparisons of these estimates with the computed probabilities for random matrices. 


III. POSSIBLE APPLICATIONS 

In the previous sections, we have showed the positive eigenpair of an irreducible symmetric non-negative matrix can 
be obtained by using an initial equal superposition state when the estimated probabilities high. This eliminates the 
necessity of knowing an initial approximate eigenvector in the applications of the phase estimation algorithm to the 
problems involving non-negative matrices. Irreducible non-negative matrices are known to have positive eigenpairs 
due to Perron-Frobenious theorem. There is also more general but the same statements indicated in this theorem for 
compact operators [11]. One may also apply permutation matrices or splitting techniques to convert the interested 
matrix to a non-negative one or to increase the expected probability outcome when the probability is low than an 
acceptable value. 

Since a wide variety of problems in science can be represented by non-negative matrices, the phase estimation 
algorithm can be applied without an initial approximate eigenvector to these problems. In the following subsections, 
we shall consider two important classes among these problems: viz., the k-local Hamiltonian problem and the stochastic 
processes. 


A. Applications to Local Hamiltonian Problems 

One of the applications of the non-negative matrices can be found in quantum mechanics where such matrices called 
stoquastic matrices [12, 13]: i.e., matrices with non-negative off-diagonal elements. If a Hamiltonian operator "H, a 
Hermitian operator acting on (C^)®", can be expressed as a sum of local Hamiltonians representing the interactions 
between at most k qubits, then it is called fc-local n-qubit Hamiltonian. More formally, TL = Finding the 

lowest eigenvalue of "H defines an eigenvalue problem which is so called local Hamiltonian problem. This eigenvalue 
problem can be also described as a decision problem to decide whether the ground state energy of H is at most a 
or at least b for given constants a and b such that a < b and b — a < l/poly{n). This problem has been shown to 
be QMA-complete for k > 2 [14] (for k=2, it is QMA-complete only when both negative and positive signs in the 
Hamiltonian are allowed [15]). 

For the local Hamiltonians with positive off-diagonal matrix elements in the standard basis, the matrix elements 
of the corresponding Gibbs density matrix, p = /Tr{—13'H), are all positive for any /3 > 0. Due to the theorem 
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given above, The ground state energy of % with positive off-diagonal elements can be associated with the eigenvector 
whose coefficients are all positive. Because one can associate a probability distribution to the ground state, the nature 
of these Hamiltonians are considered similar to stochastic processes. Thus, the Hamiltonians of these types are called 
stoquastic [12].(note that the Hamiltonian is not a stochastic matrix whose rows and columns are sum to one and the 
principal eigenvalue is one associated with an eigenvector of all ones.) 


B. Applications to Stochastic Processes 

The phase estimation algorithm can find the principal eigenpair of a stochastic matrix with probability one [9] since 
the sum of the coefficients is one for the principal eigenvector and zero for the rest of the eigenvectors. However, 
the stochastic processes are generally defined by non-Hermitian matrices, which are also widely seen in quantum 
physics: The Hamiltonian of a closed system in equilibrium must be Hermitian to have real energy values. However, 
nonequilibrium processes and open systems connected to reservoirs can be only described by non-Hermitian models 
[16]. 


1. Simulation of Non-Hermitian Operators 

Any non-Hermitian matrix "H can be decomposed into a Hermitian and a skew-Hermitian parts as: 

n = H + s=^in + n^) + ^in-n^) (i6) 

where describes the conjugate transpose of H. Here, H = -f and S = — 'H^) dehne the nearest Her¬ 

mitian and skew Hermitian matrices to %, respectively [17].Therefore, the matrix H can be used as an approximation 
to Ji in the simulation. 

A matrix TL is called normal if TL^TL — TLTL^ = 0. When "H is a normal, it is easy to see that H and S also commute: 
[H,S] = HS — HS = 0. Moreover, because of the symmetry, all the eigenvalues of H are real; and because of the 
skew-symmetry, all the eigenvalues of S involve only imaginary parts. Since TLTL^ = TL^TL, the eigenvectors of H and 
S are the same. In addition, the imaginary part of the eigenvalues of % are equal to the eigenvalues of S and the 
real parts are to the eigenvalues of H. Hence, using Ui = and U 2 = in the simulation, one can simulate the 
non-Hermitian operator "H on quantum computers by using separate two registers to obtain the imaginary and the 
real parts of the eigenvalue. However, in the stochastic cases, if TL is normal, it must be also doubly stochastic: i.e., its 
left and right principal eigenvectors are already known to be a vector of all ones with the eigenvalue one. Therefore, 
instead of an approximate normal matrix, we can use the closest Hermitian matrix defined above as H = ^{fH -T'W). 

In the non-stochastic cases, one can also try to instead find the closest normal matrix by following the Jacobi 
algorithm which attempts to diagonalize the matrix by using rotation matrices (Givens rotations), and converges to 
a matrix so called AH. The best closest matrix in the Frobenius norm is then defined by the sum of the diagonal 
elements of AH and the rotation matrices used in the algorithm [18, 19]. One can also procure a quantum circuit 
to implement this closest matrix by mapping rotation matrices to quantum gates as done in [20]. However, since 
this algorithm is based on the eigenvalue decomposition, the eigenvalues (the diagonal elements of AH, and the 
eigenvectors (the combination of the rotation matrices ) are already generated for the found normal matrix. 


IV. NUMERICAL RESULTS 

In this section, we compare the estimated af, Pregi, Preg^'- i-e. af, Pregi, Preg^: respectively, with the computed 
probabilities af , Pregi i ^reg 2 for the randomly generated matrices described below. 

A. Random Symmetric Matrices with Non-Negative Off-Diagonal Elements 

In MATLAB, we generate a random symmetric matrix with non-negative off-diagonal elements as follows: first 
a random diagonal matrix and a strictly upper triangular sparse non-negative matrix are generated. Then, these 
matrices are combined to have a symmetric matrix with non-negative off-diagonal elements: 

“/oMATLAB code to generate a random symmetric matrix 
7„ with nonnegative off-diagonal elements 
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R = triu(sprand(N,N, 0.5) , 1) ; °L N is the dimension, 

7o and 0.5 is the sparsity of the matrix. 

L = randn(N,1); 

H = R + R’ + diag(L); 

Fig.2 shows the comparisons of the estimated and computed probabilities for matrices of different dimensions generated 
by following the above code. As shown in the figure, the estimated and computed probabilities are almost the same 
and very high (close to 1). Furthermore, the success probability, Pregn with the Hadamard gates applied to \reg 2 ) in 
the end is almost one and higher than the probability, without the Hadamard gates. 


B. Random 3-Local Hamiltonians 

As our first example, the following Hamiltonian is employed: 

Pxz = ^ KijkXiXjXk + JijkZiZjZk (17) 

i,j,k 

with 0 < Kijk < 1 and — 1 < Jijk < 1- Here, X and Z are the Pauli spin operators. Choosing Kijk and Jijk randomly, 
50 numbers of random matrices are generated for each 9, 10, 11, 12, and 13 qubit systems and show the estimated 
and computed probabilities in Fig.3. 

In the second example, 1-body and 2-body interaction terms are also included separately in the Hamiltonian: 

T-Lxz = diXi + hiZi -)- KnjXiXj + JujZiZj 

* (18) 
+ K2ijkXiXjXk + J2ijkZiZjZk 
i,j,k 

with 0 < di, Knjk, K 2 ijk < 1 and —1 < hi, J^jk, J 2 ijk ^ 1- As done for Eq.(17, choosing the coefficients randomly, 
Hamiltonians of different sizes are generated. The results are shown in Fig.4. As seen in the figure, the estimation is 
not as good as the cases in Fig.2 and Fig.3. This is because the ratio defined in Eq.(15) is generally higher for the 
matrices generated from Eq.(18) in comparison to the ones from Eq.(17). 


V. CONCLUSION 

In this paper, we have showed the positive eigenpair of an irreducible symmetric non-negative matrix can be obtained 
when the estimated probabilities high. This eliminates the necessity of knowing an initial approximate eigenvector 
in the applications of the phase estimation algorithm to a wide variety of problems involving non-negative matrices. 
Here, we also have discussed how to apply it to local Hamiltonian problems and stochastic processes and showed 
the comparisons of the estimated success probabilities with the computed ones for random symmetric matrices and 
3-local Hamiltonians. 
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probabilities af, Pregi — Preg^ for random 3-local Hamiltonians of different dimensions described in Eq.(17). 
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FIG. 4. Representations of the comparisons of the estimated probabilities a?, Pregi = , Preg 2 with the computed 

Preg2 

4 

probabilities af, Pregi — Preg^ for random 3-local Hamiltonians of different dimensions described in Eq.(18). 
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